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Time structure of gamma-ray signals generated in line-of-sight 
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O ; ABSTRACT 

Blazars are expected to produce both gamma rays and cosmic rays. There- 
fore, observed high-energy gamma rays from distant blazars may contain a signif- 
icant contribution from secondary gamma rays produced along the line of sight 
by the interactions of cosmic-ray protons with background photons. Unlike the 
standard models of blazars that consider only the primary photons emitted at 
the source, models which include the cosmic-ray contribution predict that even 
~ 10 TeV photons should be detectable from distant objects with redshifts as 
high as z > 0.1. Secondary photons contribute to signals of point sources only 

if the intergalactic magnetic fields are very small, B < 10~ 14 G, and their de- 
c/3 ■ 

tection can be used to set upper bounds on magnetic fields along the line of 

sight. Secondary gamma rays have distinct spectral and temporal features. We 

explore the temporal properties of such signals using a semi-analytical formalism 

and detailed numerical simulations, which account for all the relevant processes, 



including magnetic deflections. In particular, we elucidate the interplay of time 



OO 

delays coming from the proton deflections and from the electromagnetic cascade, 
and we find that, at multi-TeV energies, secondary gamma-rays can show vari- 
es) ■ ability on timescales of years for B ~ 10~ 15 G. 
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Introduction 



Active Galactic Nuclei (AGN) are expected to be sources of both cosmic rays and gamma 
rays. While gamma rays have been observed from a number of blazars, the identification 
of cosmic rays with their sources is impossible (except at the highest energies), because the 
Galactic magnetic fields change their directions considerably However, as long as the in- 
tergalactic magnetic fields are relatively small, cosmic rays produced in blazars can travel 
close to the line of sight and produce secondary gamma rays which would significantly con- 
tribute to the radiation observed from the direction of the point sources. For nearby blazars 
such contributions are expected to be small in comparison with direct gamma-ray signals 
reaching Earth. However, for more distant blazars, the line- of-sight produced gamm a rays 
can dominate over the direct gamma rays from the source (JEssey &: Kusenkoll2010l ). The 
transition occurs because the primary gamma rays are filtered out in their interactions with 
extragalactic background light (EBL), while the fraction of secondary gamma rays produced 
by cosmic rays in in tergalactic space grows with distance traveled. Ba sed on the spectra of 
individual blazars J Essev fc Kusenkol 2010t lEssey et al.l I201CJ l2011al ) and on the trend in 
spectral softening (JEssey fc Kusenkol l201ll ). one expects the secondary contribution to be 
important for redshifts z > 0.15 and energies E > 1 TeV. 

The intrinsic gamma-ray spectra of some blazars, after correction for absorption in 
EBL, appear extremely hard, challenging the standard, e.g. the synchrotron-self-Compton 
(SSC) or External Compton models of blazars. Several solutions to this problem have been 
proposed. Intergalactic cascading of gamma rays from blazars in the case of very weak in- 
tergalactic magnetic fields (IGMFs) can in crease the effective mean free path of gamma-rays 
( lAharonian et al.ll2002t iTaylor et al.ll201ll ). however, for distant blazars this effect alone ap- 
pears to be insufficient to explain the gamma-ray spectra above 1 TeV. The very hard spectr a 
of primary gamma rays (jKatarzyhski et al.l l2006t IStecker fc Scully! 120071 ; iLefa et al.l 1201 ll ) , 
or special features in the sources ( Aharonian et al.l 120081 ) can help reconcile the data with 
theoretical pr edictions, at the cost of introducing some a d hoc assumptions. Hypothetical 
new particles (Ide Angelis et al.ll2007l ; iHorns fc Meyerll2012l ) and Lorentz invariance violation 
(IProtheroe fc Meyerl |2000| ) have been invoked to explain the data. 



The inclusion of cosmic-ray contribution offers an alternative solution to the problem. 
Indeed, since a significant fraction of gamma rays in this model is produced relatively close 
to the observer, this model reduces dramatically the impact of gamma-ray absorption in 
EBL. This is illustrated in Fig. [TJ Protons with energy E > 10 16 eV propagating through 
weak IGMFs without strong deviations from the line of sight can carry energy from the 
source close to the observer and can generate a substantial gamma-ray flux at multi-TeV 
energies. Remarkably, the predicted spectra of secondary gamma rays depend only on the 
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source redshift (determined from independent observations). For each source, the power 
emitted in cosmic rays is the only fitting parameter which can be used to fit the data. As 
long as the source redshifts are known, the predictions of this model are solid and robust. 
The spectra calculate d for the redshifts of all observed distant blazars provide a very good fi t 
to observational data (jEssey fc Kusenkoll2010l ; lEssey et al.ll2010l l2011al ; iMurase et al.ll201ll ). 
with a reasonable required luminosities in cosmic rays a ssuming that the escape of protons 



from the source is strongly beamed toward the observer (jEssey et al.ll2011al ; iRazzaque et al. 
20121 ). 



Confirmation of this model by future observations will have several important conse- 
quences. It will imply that (i) cosmic rays are, indeed, accelerated in AGN, as has long been 
suspected, but never before confirmed by observations; (ii) interga lactic magnetic field s are 
fairly small, of the order of several femtogauss (10~ 15 G) or lower (jEssey et al.ll2011bl ); and 
(iii) the problem of intergalactic gamma-rays can be somewhat relaxed, and consequently 
the upper limits on EBL derived while neglecting the cosmic-ray contributions may need be 
revised. Within this model, the expected temporal structure of signals from distant blazars 
at the highest energies should reflect the time delays cosmic rays undergo in the intergalac- 
tic magnetic fields. We note that while time variability has been observed for nearby TeV 
blazars at TeV energies and for distant TeV blazars at energies above a few hundred GeV, 
no variability has been reported so far for distant TeV blazars at TeV energies. Here we 
call distant those blazars that have large enough redshifts for the primary TeV gamma-rays 
to die out. Based on the spectral fits (jEssey fc Kusenkol l2010t Essev et al.l |2010| . l2011al : 



Murase et al.ll201l[ ). and the spectral softening of blazars lEssey fc Kusenkol (J20111 ). one con- 
cludes that most blazars with redshift z > 0.15 should fall in this category. Since the ratio of 
gamma-ray luminosity to cosmic-ray luminosity can vary from source to source, one expects 
a mixed population to exist at some intermediate redshifts 0.1 < z < 0.15, where stronger 
cosmic ray emitters should be observed in secondary gamma rays, while stronger gamma ray 
emitters should be observed in primary gamma rays. Furthermore, if stronger IGMFs exist 
in the directions of specific sources, the secondary contribution can be suppressed. For ex- 
ample, PKS 2155-304 at z — 0.12 is an example of a source at an inter mediate redshift from 



which primary signals are observed, as indicated by its TeV variability IHESS Collaboration 



( j2010l ). Whether the lack of TeV variability is a generic feature of distant blazars, or merely 
an artifact of low statistics in multi-TeV photons, should be clarified in future observations. 
An important issue in this context is the knowledge of the spectral and temporal features 
of the radiation pre dicted by the mode l. The spectral features of the radiation have been 
studied in detail by lEssey et al.l ( j2011af h 



In this paper we consider the extent to which the time variability at high energies should 
be erased by the cosmic ray propagation delays. We will focus on calculating the Green's 
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Fig. 1. — Secondary gamma rays produced in line-of-sight interactions of cosmic rays result 
in harder spectra for distant sources. Since most of the observed photons are produced 
relatively close to the observer, there is less attenuation due to the interactions with EBL. 
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function, which corresponds to a time delay from an infinitely narrow pulse of protons at 
the source. Realistic time profiles can be obtained by convolving the time-dependent source 
luminosity with this Green's fucntion. However, since a fair fraction of the blazar flaring 
activity occurs on the time scales much sorter than those we discuss below, in many cases the 
Green's function can be interpreted as a distribution of photon arrival times from a flaring 
source. In applications of our method to data analysis, one can employ time-dependent 
templates inferred from lower energies. 



2. Basic estimates and scaling laws 

There is no doubt that, for nearby blazars, the primary gamma rays produced at the 
source are responsible for most of the observed radiation. While these objects can also 
produce cosmic rays, the contribution of secondary gamma rays is not expected to dominate. 
However, for larger distances, the primary gamma-ray component is filtered out above a TeV, 
while the secondary contribution is enhanced. Indeed, the scaling of the primary gamma rays 
with distance is determined by the losses due to pair production in gamma-ray interactions 
with extragalactic background light (EBL): 

-Fprimary^eO OC — eXp(-d/A 7 ). (1) 

In contrast, gamma rays gener ated in line-of-sight interactions of cosmic rays ex hibit a very 



different scaling with distance (jEssey fc Kusenkoll2010l : lEssey et al.ll2010l l2011ai ): 



A., 



secondary 



r \ 

(2) 



1/d, for d <C A 7 , 
1/d 2 , for d > A 7 . 

Here A 7 is the distance at which EBL opacity to TeV gamma rays is of the order of 1. The 
lack of suppression is due to the fact that the photon backgrounds (CMB and EBL) act as 
a target on which gamma rays are produced by the cosmic rays. Hence, a higher column 
density of background photons for a more distant source boosts, not hinders the gamma-ray 
production. 

As long as IGMFs are weak enough to cause only small deflections, for a sufficiently 
distant source, secondary gamma rays dominate because they don't suffer from exponential 
suppression as in Eq. (JTJ), which is absent from Eq. (J2J). The transition from primary to sec- 
ondary photons occurs when the optical depth to pair production exceeds 1 . The correspond- 



ing redshift can also be inferred from the spectral softening of the spectra (JEssey fc Kusenko 
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201ll ). Based on these estimates, one can expect the secondary gamma rays to dominate for 
sources at redshifts z > 0.15. 



The success of the s pectra l fits to the data for secondary gamma rays (JEssey fc Kusenko 



2010l : lEssey et al.ll2010l . l2011ai ) can be interpreted as possible evidence of cosmic ray accel- 
eration in blazars. Within this interpretation, the beamed energy output in E > 10 17 eV 
cosmic rays required to fit the o bserved spectra of distant blazars is of the order o f 10 43 erg, 



or 10 45 erg isotropic equivalent (lEssev fc Kusenkoll2010l ; lEssey et al.ll2010l . l2011al ). which is 
consistent with many models ( IBerezinsky et al.l 120061 ) . The luminosity required to explain 
ultrahigh-energy cosmic rays (UHECR) depends on the assumed spectrum, which is often 
parameterize d by a broken power law with a break at some value E c (unknown a priori). 
According to IBerezinsky et al.l (120061 ) . the AGN luminosities in cosmic-ray protons needed 
to account for the UHECR data are 5.6 x 10 43 erg/s, 2.5 x 10 44 erg/s, and 1.1 x 10 45 erg/s, for 
E c = 10 18 eV, E c = 10 17 eV, and E c = 10 16 eV, respectively. These estimates are in general 
agreement with our results. We also note that, because of the selection effects, the sources 
observed from large distances are by no means average: they are the brightest AGNs, which 
generate exceptional power in cosmic rays. 

The spectra of observed gamma rays generated in this fashion depend on the intervening 
intergalactic magnetic fields. It is easy to understand some qualitative features of this 
dependence in terms of a simplified random- walk description of the proton propagation. Let 
us consider a short pulse of protons emitted from a source at distance d. At later times, the 
pr oton pulse broadens and takes the shape f(t, r). The explicit fo rm of f(t, r) was compute d 
bvlAharonian et al.l ( ]2010l ) and will be discussed below. (See also lAlcock fc Hatchettl ( 119781 ); 
Williamson! Jl972l 



At every point in its trajectory, the proton interactions with the cosmic background 
generate a flux of gamma rays, which quickly (on a kpc length scale) cascade down to 
energies below the threshold. From that point on, gamma rays travel without further time 
delays. However, during the cascade development, the IGMFs cause some delays (which are 
longer than the delays of the protons in the IGMFs for E^ below 10 TeV). 

Let us consider the proton propagation in IGMFs. We assume that IGMFs form a 
lattice with correlation length l c , in which a proton with energy E p = 10 17 eV random-walks 
over a distance d ~ 1000 Mpc = n x l c , where n ~ 10 3 . 

The angle between the proton momentum and the line of sight performs a two-dimensional 
random walk in small steps of ~ 10~ 6 . When the proton interacts with a background photon, 
it emits a narrow shower in the direction of the proton's momentum. The prompt gamma 
rays are emitted into a narrow angle (~ -E p /MeV) _1 for Bethe-Heitler pair production, or 
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(~ Ep/0.2 GeV) -1 for pion photoproduction. The cascade develops and broadens this angle, 
with larger angles at lower energies. To shower in the direction of observer after n steps of 
random walk, the proton angle should return to 0. For a 2D random walk, the probability of 
not returning to zero after n steps is 72(71) = ^ + O ( (ln \ 2 j . This probability drops below 
1/2 for n > exp(27r) ~ 5 x 10 2 . For d ~ 1000 Mpc, n ~ 10 3 , and so each proton angle returns 
to the origin about ~ 1 time per distance traveled. Therefore, the diffusion approximation 
is justified, and a "typical" delay can be computed using the distance traveled, assuming the 
random walk. 

Deflection of a proton in a single cell is 6q. This deflection and the time delay are 
determined by the Larmor radius 

R B = -^ = 10 5 Mpc E pA7 /B_ 15 , 

where E Pj n is the proton energy in units of 10 17 eV, and -8-15 is the value of the magnetic 
field in femtogauss. 

Therefore, 



* = £ = 10 " (wc ' B - IE > 
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Time delay in crossing a single cell is 



At = hi = 10 4 s ( B - 15 



c u \E 



After n ~ 10 3 steps of random walk, the time delay is 

/ ' B \ 2 

r p = nAt ~ 10 7 s — - . (3) 

Significant time delays are also incurred in the EM showering process. Each observed 
gamma ray was at some point an electron in the cascade. The time delay of each gamma ray 
with an observed energy E 1 is dominated by the delay during the lowest-energy "electron" 
phase of this gamma ray. The electron energy is related to the energy of observed IC 7 ray 
by Ky oc E\. IGMFs act on the electron over a distance of the order of its cool ing distance 



D P oc 1/E P . The time del ay incurred in this process is proportionate to the sum (jlchiki et al. 



20081 ; iMurase et al.l 120081 ) of D e and the mean free path to pair production App, which has 
the same energy dependence (with a much larger prefactor), App oc 1/E e . The resulting 

delay is 

eB 
T e = (App + D e )9l/c, where 9 = D e —. 



Therefore, 



n2n ^,e 2 B 2 B 2 B 2 

r y ^r e = D 2 e (X PP + D e )—— oc — oc —— (4) 



cE 2 E B E ,/2 

where we have assumed D e ~ kpc <C l c . 



The proton delay ([2]) is 



r p oc — . (5) 

p 



The total time delay of an observed photon is the sum of r p and r 7 : 



^ 2 E 5 J r 



r tot — r p + r 7 — Ci^r + C2 „ K / 9 5 (6) 



where Ci and C2 are some constants. 



One can, therefore, expect the following structure of time delays. The shortest delay 
time is determined by the delay in the arrival of the highest-energy proton; this time delay 
is given by Eq. fl5]). High-energy protons travel faster than the gamma-ray cascades, and 
they are followed by a tail of trailing gamma rays. There are two contributions to the total 
delay time, which have a different dependence on energy (^. As the second term in Eq. §§§ 
diminishes with energy, the time delay approaches a plateau independent of the photon 
energy. The height of this plateau is determined by the energy of the proton. This agrees 
with the results of numerical caclculations presented in Figure |2j 



3. Semi-analytical description 

Let us consider the time delays due to the propagation of protons. Protons interactions 
with extragalactic background radiation (both CMB and EBL) occur with a very low prob- 
ability for energies below the pion production threshold. In the pair production process, a 
proton loses only ~ 10~ 3 of its energy in each collision. Thus, one can neglect the energy 
losses for protons in making some basic estimates. (However, in our numerical calculations, 
we take into account all the energy losses, including adiabatic losses.) Also, the deflec- 
tion angles can be assumed small for the relevant range of parameters. Calculation of the 
electromagnetic cascade initiated by the secondary gamma rays produced in proton-photon 
interactions is much more difficult, and there is no simple analytical approach that could 
allow one to calculate the distribution function of gamma rays. Therefore we will employ 
a hybrid approach by combining an analytical treatment of protons with a Monte Carlo 
simulation for the electromagnetic cascade. 
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To calculate the distribution function of protons, let us consider a mono-energetic beam 
of protons emitted with energy E at some point in time. Random deflections in weak 
IGMFs result in arrival time distribution which is convenient to consider as a function of 
time delay parameter r = t — r/c, where r is the distance to the source and c is the speed of 
l ight. In a small- angle approximation, one can express the distribution function as follows 



(jAharonian et al.l 120101 ): 



I / CT ~ J CT \ \ 

UE ^ r) = -r\-FW) fA \^W)))' (7) 

where 

oo 

h(v) = 4vr 2 ^(-lrVr^ (8) 



'|2\ 



n=l 

l 



;> = |(|) a?). ( 9) 

Here (8 2 ) is the mean square deflection angle per unit length. The correlation length l c and 
the mean square of magnetic field (B 2 ) enter in Eq. (J7J) as parameters. The normalization 

oo 

of the function Ja is so that J JaAt = 1- Then the distribution function of protons injected 

o 
with energy spectrum J P {E) at the distance r from the source is 

fr[EiTiT)= Msm^i. (10) 

Let protons interact with low energy photon field f p h{e)- The protons with a monoen- 
ergetic distribution (normalized to one particle) produce electron-positron pai rs at the rate 



$(EL E p ), where E p is the energy of protons, E e is the energy of pairs. Following lKelner fc Aharonian 
(120081 ). one can express Q(E e ,E p ) as follows: 



/rjC) k ■ II r 

d^f P h(e)^J 8{E e -c(u irPe ))da, (11) 

where k, u p , and uif are four- velocities of photon, proton, and the laboratory frame, respec- 
tively; and p e is four- momentum of electron (or positron). 7 e and e are the proton Lorentz 
factor and the energy of photon in the laboratory frame, da is the Bethe-Heitler cross sec- 
tion. The photon field used in this calculations includes the CMB and EBL, and one can 
neglect the redshift evolution. 

Then distribution function of electrons produced at the distance r with inherited time 
delay r is 

/.(£., r,r) = fdE / p{Ep)f f 2 Ep ' T ' r) <l>(E e ,E p ). (12) 
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These electrons initiate an electromagnetic cascade. Let f cas {E e ,E^,s) be the number of 
photons with energy E 1 produced in the cascade initiated by an electron with energy E e at 
the distance s from the observer and detected at the point of observation. The mean time 
delay of the photons is r cas (E e , E 7 , s). The extragalactic magnetic field is a parameter. Then 
for the UHE proton source at the distance d the number of photons produced in the cascade 
with time delay r = r cas + T prot is 

d 

f 7 (Ej,r,d) = dr dE e f e (E e ,T -T cas (E e ,E^,r),d-r)f cas (E e ,E 1} r) (13) 

o 

The Eqs. ([7|), f fl2|) and ffT3|) constitute the integral which is calculated numerically. The 
functions f cas and r cas actually depend on redshift but not on distance, therefore we use the 
relation 

dr = 4-r ^=dz, (14) 

#o (i + z)y/(i + z ) 3 n m + n A 

to express distance via redshift and perform integration over z. 

We assume that the source produces a power-law spectrum of protons with a spectral 
index a and the energy range from 0.1 Eq to E$. The results for the mean time delay of gamma 
rays are presented in Figs. |2] and [3] as a function of the varying cut-off energy, distance to the 
source and spectral index. In these calculations we assumed that the intergalactic magnetic 
field has the strength B = 10~ 15 G and the coherence length l c — 1 Mpc. Unless specified 
otherwise, we use E = 10 17 eV, a = 2. 

The advantage of the analytical description presented above is the possibility to study 
the time delay distribution of gamma rays for a variety of initial proton spectrum parameters. 
The numerical Monte Carlo approach described below has computational limitations on the 
number of initial particles, which can complicate the study of how a proton spectrum with 
a wide energy range can affect the time delays of gamma rays. The protons injected with 
slightly different energies would have the time distribution which is similar to the time 
distribution for a monoenergetic proton beam. In contrast, the time distribution of protons 
with a broad energy spectrum is a sum of time distributions of protons with different energies, 
which is stretched out in time. This would spread the arrival times of gamma rays along 
a large time span. The illustration of this effect can be seen from the comparison of the 
second panel of Fig. H] and Fig. EJ The range of the proton energies has a strong effect on 
the lower energy gamma rays, whereas the time distribution for E^ — 1 • 10 14 eV does not 
change significantly. This is because, in the case of a broad spectrum, protons with different 
energies can contribute gamma rays of a given energy. On the other hand, only the protons 
of highest energies are responsible for the production of gamma rays of E~ = 1 • 10 14 eV. 
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Fig. 2. — Mean time delay of gamma rays at redshift z = 0.17 for different cutoff energies 
E of proton spectrum. 
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Therefore, the corresponding time distribution is similar to the one for the monoenergetic 
protons (c/. Fig. [6]). 

The flux of gamma rays arriving at any given time comprises contributions from protons 
at different points along the line of sight. We obtain this flux by integrating over the proton 
distributions shifted by a time delay incurred in the electromagnetic cascade. The latter 
delays were obtained from numerical calculations using a single Monte Carlo numerical run. 
The results are shown in Figs. S] and Fig. The multiple-peak structure apparent in these 
curves is the result of adding contributions from different distances with the delay profile 
obtained from a single numerical run. If we adopted a different approach and used an 
averaged delay profile, as in Fig. [7J the "many-peak" structure would be erased. 



4. Numerical Monte-Carlo calculations 

In addition to the semi-analytical results we also performed a full scale Monte Carlo 
simulation to track the arrival times of individual particles. The source was modeled by 
an instantaneous pulse of protons to represent the Green's function needed to calculate 
the distribution of arrival times. Particles are advanced in time steps of roughly 0.1 — 
1 kpc, updating momentum, position and time delay with all relevant interactions taken 
into account. Gamma rays arriving at the z = surface are binned and the mean arrival 
time and standard deviation are calculated. 



The proton energy loss processes are well studied (ISzabo fc Protherodll994j ) and can 
be described by a standard approach. We calculate all the relevant energy losses, includ- 
ing adiabatic losses and the losses due to the interactions with photon backgrounds. The 
most important contributions to secondary photon production are photopion production and 
proton pair production (PPP). 

The photopion production processes involve the following reactions: 

V + lb ->■ Tl + 7T + 

P + lb ->■ P + TT° (15) 

where 7^ is either a CMB or EBL photon. PPP occurs in the reaction 

P + lb ~^p + e + + e~. (16) 

The pair production on the CMB is the dominant reaction, but pion photoproduction on 
EBL also contributes. Pion photoproduction on CMB has a threshold above 10 19 eV, but 
pion production on EBL is possible for all energies we consider. The efficiency of energy 
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transfer to the electromagnetic shower depends on the proton energy and on the dista nce to 
the source. A more detailed discussion is presented elsewhere (jAharonian et al.ll2012[ ). 



The mean interaction length, A, for a proton of energy E traveling through a photon 
field is given by 



[A^]]- 1 



8(3E 2 



00 n(e^ r s ™**( e > E ) 



cr(s)(s — mfydsde, 



(17) 



where n(e) is the differential photon number density of photons of energy e, and a(s) is the 
appropriate total cross section for the given process for the center of momentum (CM) frame 
energy squared, s, given by 



s = ml + 2eE(l - (3 cost 



'11 



where 9 is the angle between the proton and photon, and (3 is the proton's velocity. 
For pion photoproduction, 

s min = (ml + ml) 2 (19) 



and 

For proton pair production 
and 

For both processes, 



m^im^ + 2m p 
2E(1 + P) 



Smin = \m v + 2m 



2\2 



m e (m e + m p ) 

E 



s max (e,E) = ml + 2eE(l + (3). 



(20) 
(21) 
(22) 
(23) 



Both pions and neutrons quickly decay via the processes 



n — > p + e~ 



v, 



IX 



IX 



-+ /i + + Vp -> e + + v e + Vp + u^, 

->■ 2 7 . 



(24) 



The outgoing dis tribution functions for pion photoproduction were generated using the 
SOPHIA package iMiicke et all J200nh . 



Primary gamma rays and gamma rays produced from the above equations can interact 
and pair-produce on background photons. The resulting electron positron pairs will IC 
scatter CMB photons. The upscattered photons can once again pair produce, this chain 
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reaction is known as electromagnetic (EM) showering. 

The interaction length for photons for pair production off the EBL is 



[A]" 1 = (^f ) J^ e" 2 n(e) J ™ e 2sa(s)dsde, (25) 



where 



and 



(^) 2 (1 - /3 2 ) [(3 - /3 4 ) In 1±| - 2/3(2 - /3 2 )] (26) 



1 f e2 Yn 

cr = -n{ — ) (I 



/3 = (1 - 1/,) 1/2 , (27) 

and n(e) is the differential photon number density of photons of energy e. 

To simulate magnetic field effects, the IGMF is modeled by cubic cells of a given mag- 
netic field strength with sides equal to the chosen correlation length, l c , and a random 
direction. Particles are moved forward in fine time steps and the deflection of the particle is 
calculated using the Larmor radius and IGMF direction. Time delays for charged particles 
are calculated in comparison to a photon traveling in a straight line to the observer. 

For the analysis of time delays, we have performed multiple runs and averaged the 
results, as shown in Fig. 

The results of the simulation are shown in Fig. [BJ where delays from deflections in the 
IGMF are shown for a source at z = 0.2. It is evident that secondary photons produced at 
large distances conform to the power-law behavior as in Eq. (J41) . This approximate power law 
is illustrated in Fig. [H] by a dashed line. The flattening at low energies can be understood by 
the way the code handles deflections. Particles are moved forward in time steps of roughly 
0.1 — 1 kpc and deflections are assumed to be less than tx within a single time step. For the 
lowest energies this is not always true and the code will underestimate the deflection, thus 
producing time delays below the power law. 



5. Discussion and conclusion 

The main qualitative features of the Green's function computed numerically and shown 
in Fig. [2] and Fig. [HJcan be easily understood. For lower energies (below TeV), time delays 
r oc B 2 E~ 5 / 2 d, where d is distance to the source. The nearby showers arrive before distant 
showers, so that the late arriving gamma rays have lower energies and longer delays. The 
plateau that develops at E > ITeV is due to the prompt showers emitted by the protons 
nearby, for which the time delays are determined by the proton deflections in IGMFs. In the 
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absence of cosmic rays, the spectrum would drop above 1 TeV, and the multi-TeV gamma 
rays would not be observed. 

For sub-TeV secondary gamma rays the electromagnetic cascade delays are always longer 
than the proton delays, and the arrival photons peak at the time given by Eq. (J3J). At energies 
above TeV, the proton delays come to dominate, in accordance with the broken power law 
in Eq. (|6]). The numerical results differ somewhat from the scaling in Eq. ([6]). In particular, 
at low energies, the delays appear to scale as E~ 2 rather than E~ 25 . The difference can be 
explained by a combination of several effects. For electron energies below 30 GeV, the cooling 
distance exceeds the magnetic field correlation length, which we assumed to be l c ~ 1 Mpc. 
This changes the energy dependence in Eq. (TjJ because the energy-dependent cooling distance 
must be replaced by the constant correlation length. Furthermore, integration over energies 
in the cascade affects the power-law behavior. For these reasons, our basic estimates in 
section 2 were not expected to capture all the features evident in the numerical results. 

The proton delay is strongly dependent on the high energy cutoff of the cosmic ray 
source, which affects the energy at which the proton delays begin to dominate. This can be 
seen in Fig. |HJ This behavior is further illustrated in Fig. H] and Fig. EJ where one can see 
that the proton delays begin to dominate at E ~ 10 TeV for a proton high energy cutoff of 
10 8 GeV and an IGMF = 10~ 15 G. 

The distribution of gamma-ray arrival times depends on the injection spectrum of pro- 
tons, as one can see from a comparison of Fig. H] and Fig. |5J This is in contr ast with the 



spectra of ga mma rays, which are not sensitive to the proton injection spectrum (lEssey et al. 
2010l . l2011al ). Hence, one can, at least in principle, learn about the proton injection spectrum 
from timing observations, but not from the spectra alone. (Neutrino spectra also depend 



on the proton injection spectrum ( lEssey et al.ll2010l ).) Furthermore, stochastic broadening 



" 



illustrated in Fig. [7] also affects the predictions. 

Based on our results, the observed time variability should be washed out on time scales 
shorter than ~ 0.1 yr, for distant blazars (z > 0.2), at TeV and higher energies. Time 
variability can be present for z > 0.2, E > 1 TeV on the time scales of 0.1 — 1 yr. If gamma 
rays with E ~ 10 2 TeV are observed, they can exhibit variability on shorter time scales. 

Of course, one must also consider the delays the cosmic rays undergo at the source. 
Blazars are known to be highly variable, and this variability could affect the shape of the 
observed spectrum. The magnetic fields within galaxies are on the order of 1 /iG which can 
lead to significant delays in the source. On the other hand, the structure of magnetic fields 
in front of the blazar jets is not known. Furthermore, the effect of the source variability 
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would be to suppress the observed power of the source by a factor 

/damp ~ ^active ( ~ ] , (28) 

\ Mclay / 

where tdeiay is the typical proton delay at the source, t ac tive is the typical time the source is 
active or flaring and iV active is the number of times th e source is active in the time period 



tdeiay This damping should not be a significant effect (jDermer et al.ll201ll ). especially since 
the typical deflections at the source are not big enough to affect the beaming factors assumed 
in popular models. 

An alternative situation is that the magnetic fields within the blazar jet are not randomly 
distributed, but are, instead, strongly correlated with the direction of the jet. Blazar jets 
emit an extremely large amount of charged matter and the wind in the direction of the jet 
can eliminate any random-field configuration that one usually expects in a galaxy. Thus, 
it is possible that cosmic rays escape the source along the jet with very small time delays, 
preserving the intrinsic variability of the source. In this case, delays in the intergalactic 
medium can broaden the intrinsic variability to the energy dependent timescales of these 
delays. 

At energies where the optical depth of the observed gamma rays is below one, we expect 
the signal to be primary gamma rays, and any variability in this signal is indicative of the 
source variability. This variability should not depend strongly on the energies of the gamma 
rays, but rather on the scale of the structure at the source producing the gamma rays. 

However, we expect a very different behavior for energies at which primary gamma rays 
are significantly attenuated by pair production off EBL. In the case of strongly correlated 
magnetic fields in the jet, we expect that the variability should show different structure in the 
low energy component, where it should depend on the energy. The spectrum should show 
variability on shorter timescales for higher energies until some critical energy E c ~ TeV, 
where the timescales cease to decrease further, thanks to the domination of the cosmic ray 
contribution. In the case of large delays within the source, we expect all variability to be 
washed out at these higher energies (typically around a TeV for most observed sources). 

Some exceptionally bright flares can come through around E c and rise above the pedestal 
created by the stochastic arrival times of protons. Such flares should have distinctly softer 
spectra than the hard pedestal, which can be a means of distinguishing these flares from the 
stochastic pedestal. 

For most of discussion, we have assumed that IGMFs have strengt hs are of the order o f 



a femtogauss. This range is suggested by the spectral fits to the data lEssey et al.1 (j2011bl ). 



However, field strengths well below a femtogauss can be consistent with the data as well. 



-17- 



In the case of very weak IGMFs, the time delays become smaller, since r oc B 2 . For 
B ~ 10~ 18 G, the time delays can be as short as minutes. 

We have also assumed that the strength of IGMF is constant on average, which is a good 
assumption for propagation in the voids. However, if the line of sight intersects a filament 
of stronger, e.g., nanogauss magnetic field, the reduction of the secondary photon signal 
depends on the size of the filament and on its location. Thin filaments can only intercept as 
small fraction of protons within the 0.1 degree associated with a given source. However, a 
thick filament or a sheet of strong field can deflect protons, reducing the secondary signal. 

Temporal structure of gamma-ray signals can be used to measure the IGMF structure 
and EBL intensity in different directions, on a source-by-source basis. In addition, it may 
provide a way to probe the high energy cutoff of cosmic ray sources, as well as the spectrum 
of EBL. A statistical analysis on multiple bins of data is needed to determine the variability 
at different energies. This presents challenges at the highest energies because of the low 
statistics, but longer observation times and the advent of next generation experiments should 
make this analysis increasingly powerful. 
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Fig. 3. — Left panel: mean time delay of gamma rays for the sources at different redshifts 
and the proton spectrum with cutoff energy Eq = 10 17 eV. Right panel: mean time delay of 
gamma rays for the source at z = 0.17 and the proton spectrum with cutoff energy Eq = 10 17 
eV and different spectral indices a. 
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Fig. 4. — Time delay distribution of gamma rays in arbitrary units (the maximum of dis- 
tribution is normalized to unity) at different energies from the source at z = 0.17. Each 
plot corresponds to the proton spectrum with different cutoff energy E$ and spectral index 
a = 2. The injected spectra of protons are taken in the range from O.IEq to E . 
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Fig. 5. — Time delay distribution of gamma-ray in arbitrary units (the maximum of distri- 
bution is normalized to unity) from the source at z — 0.17. The injected spectrum of protons 
is almost monoenergetic with energy E = 10 17 eV. 
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Fig. 6. — Arrival time probability distribution in arbitrary units for secondary gamma rays 
with energies 1 TeV (blue, long-dashed line), 10 TeV (purple, short-dashed line) and 100 TeV 
(red, solid line). Results are shown for a cosmic ray source at z — 0.2 with a high energy 
cutoff of 10 8 GeV and an IGMF of 10" 15 G. 
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Fig. 7. — Arrival time probability distribution in arbitrary units for 1 TeV secondary photons 
for multiple numerical runs. The results shown are for roughly 300,000 secondary photons 
with an IGMF of 10~ 15 G and UHECR cutoff of 10 10 GeV. The blue thick line represents 
the sum of all distributions and the thin red lines are a representative set of distributions. 
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Fig. 8. — Arrival time delays from an instantaneous pulse emitted by a source at z — 0.2, 
assuming B = 10~ 15 G with l c — 1 Mpc correlation length. 
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Fig. 9.— Arrival time probability distribution in arbitrary units for primary cosmic rays. 
Results are shown for a cosmic ray source at z = 0.2 with a high energy cutoff of 10 10 GeV 
(blue, long dashed), 10 9 GeV (purple, short dashed) and 10 8 GeV (red, solid) and an IGMF 
of 1(T 15 G. 



